Example Solow model
Contents
3. Example Solow model¶
In this jupyter notebook we will specify, solve and analyse a simple Solow model in ModelFlow.
#Required packages
import pandas as pd
# Modelflow modules
from modelclass import model
#for publication
latex=0
model.widescreen()
3.1. Specify the model¶
We start by defining the logic of the Solow model in the Business Logic Language.
fsolow = '''\
Income = a * Capital**alfa * Labor **(1-alfa)
Consumption = (1-saving_rate) * Income
Investment = Income - Consumption
diff(Capital) = Investment-Depreciation_rate * Capital(-1)
diff(Labor) = Labor_growth * Labor(-1)
Capital_intensity = Capital/Labor
'''
3.2. Create a model class instance¶
After defining the Business Logic Language and storing it in the variable ‘fsolow’, we create a class instance called msolow.
msolow = model.from_eq(fsolow,modelname='Solow model')
Above the equation for Capital and Labor on the left hand side of the = (equal to) consist of an expressing diff(Capital) and diff(Labor). The equations are not normalized.
To solve a model in modelflow all equations has to be normalized. Meaning that the left hand side only consist of variables not expressions. So the function model.from_eq will normalize the equations as the first step before the model can be solved.
In this case first diff(Capital) is transformed to \(\Delta capital = capital-capital(-1)\). Then the lagged variables is moved to the right side of the =.
The same goes for diff(labor).
So the normalized business language of the model now looks like:
msolow.print_model
FRML <> INCOME = A * CAPITAL**ALFA * LABOR **(1-ALFA) $
FRML <> CONSUMPTION = (1-SAVING_RATE) * INCOME $
FRML <> INVESTMENT = INCOME - CONSUMPTION $
FRML <> CAPITAL=CAPITAL(-1)+(INVESTMENT-DEPRECIATION_RATE * CAPITAL(-1))$
FRML <> LABOR=LABOR(-1)+(LABOR_GROWTH * LABOR(-1))$
FRML <> CAPITAL_INTENSITY = CAPITAL/LABOR $
3.3. Create some data¶
To show what Modelflow can do, we create a Pandas dataframe with input data. And print the first 5 out of 300 observations.
N = 300
df = pd.DataFrame({'LABOR':[100]*N,
'CAPITAL':[100]*N,
'ALFA':[0.5]*N,
'A': [1]*N,
'DEPRECIATION_RATE': [0.05]*N,
'LABOR_GROWTH': [0.01]*N,
'SAVING_RATE':[0.05]*N})
df.head(2) #this prints out the first 5 rows of the dataframe
| LABOR | CAPITAL | ALFA | A | DEPRECIATION_RATE | LABOR_GROWTH | SAVING_RATE | |
|---|---|---|---|---|---|---|---|
| 0 | 100 | 100 | 0.5 | 1 | 0.05 | 0.01 | 0.05 |
| 1 | 100 | 100 | 0.5 | 1 | 0.05 | 0.01 | 0.05 |
3.4. Run the model¶
result = msolow(df,keep='Baseline') # The model is simulated for all years possible
result.head(29)
| LABOR | CAPITAL | ALFA | A | DEPRECIATION_RATE | LABOR_GROWTH | SAVING_RATE | INVESTMENT | CAPITAL_INTENSITY | INCOME | CONSUMPTION | |
|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 100.000000 | 100.000000 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 0.000000 | 0.000000 | 0.000000 | 0.000000 |
| 1 | 101.000000 | 100.025580 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.025580 | 0.990352 | 100.511609 | 95.486029 |
| 2 | 102.010000 | 100.076226 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.051924 | 0.981043 | 101.038487 | 95.986562 |
| 3 | 103.030100 | 100.151443 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.079029 | 0.972060 | 101.580575 | 96.501546 |
| 4 | 104.060401 | 100.250762 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.106891 | 0.963390 | 102.137821 | 97.030930 |
| 5 | 105.101005 | 100.373733 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.135509 | 0.955022 | 102.710176 | 97.574667 |
| 6 | 106.152015 | 100.519926 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.164880 | 0.946943 | 103.297593 | 98.132713 |
| 7 | 107.213535 | 100.688931 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.195002 | 0.939144 | 103.900030 | 98.705029 |
| 8 | 108.285671 | 100.880357 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.225872 | 0.931613 | 104.517449 | 99.291576 |
| 9 | 109.368527 | 101.093830 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.257491 | 0.924341 | 105.149813 | 99.892323 |
| 10 | 110.462213 | 101.328993 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.289855 | 0.917318 | 105.797092 | 100.507238 |
| 11 | 111.566835 | 101.585506 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.322963 | 0.910535 | 106.459257 | 101.136294 |
| 12 | 112.682503 | 101.863045 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.356814 | 0.903983 | 107.136282 | 101.779468 |
| 13 | 113.809328 | 102.161300 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.391407 | 0.897653 | 107.828145 | 102.436738 |
| 14 | 114.947421 | 102.479976 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.426741 | 0.891538 | 108.534829 | 103.108087 |
| 15 | 116.096896 | 102.818793 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.462816 | 0.885629 | 109.256317 | 103.793501 |
| 16 | 117.257864 | 103.177483 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.499630 | 0.879920 | 109.992597 | 104.492967 |
| 17 | 118.430443 | 103.555792 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.537183 | 0.874402 | 110.743661 | 105.206478 |
| 18 | 119.614748 | 103.953478 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.575475 | 0.869069 | 111.509502 | 105.934027 |
| 19 | 120.810895 | 104.370310 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.614506 | 0.863915 | 112.290118 | 106.675612 |
| 20 | 122.019004 | 104.806070 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.654275 | 0.858932 | 113.085509 | 107.431233 |
| 21 | 123.239194 | 105.260550 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.694784 | 0.854116 | 113.895678 | 108.200894 |
| 22 | 124.471586 | 105.733554 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.736032 | 0.849459 | 114.720631 | 108.984599 |
| 23 | 125.716302 | 106.224895 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.778019 | 0.844957 | 115.560378 | 109.782359 |
| 24 | 126.973465 | 106.734397 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.820747 | 0.840604 | 116.414931 | 110.594185 |
| 25 | 128.243200 | 107.261893 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.864215 | 0.836394 | 117.284305 | 111.420090 |
| 26 | 129.525631 | 107.807224 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.908426 | 0.832323 | 118.168518 | 112.260093 |
| 27 | 130.820888 | 108.370242 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.953380 | 0.828386 | 119.067591 | 113.114212 |
| 28 | 132.129097 | 108.950808 | 0.5 | 1.0 | 0.05 | 0.01 | 0.05 | 5.999077 | 0.824578 | 119.981548 | 113.982470 |
3.5. Create a scenario and run again¶
dfscenario = df.upd('LABOR_GROWTH + 0.002') # create a new dataframe, increase LABOR_GROWTH by 0.002
scenario = msolow(dfscenario,keep='Higher labor growth ') # simulate the model
3.6. Now the results are also embedded in msolow.¶
.basedfcontains the first run of the model.lastdfcontains the last run of the model
Also in this case the keyword keep is used. This causes the results to be stored in a dictionary msolow.keep_solutions. This can be useful when comparing several scenarios.
3.7. Inspect results¶
3.7.1. Using the [ ] operator¶
We can select the variables of interest with wildcards. This will operate the results stored in basedf and .lastdf
3.7.1.1. Look at variables starting with a C¶
msolow['#ENDO']
3.7.1.2. Look at all endogenous variables¶
msolow['labor*'].dif.plot()
3.7.2. Using the keept solutions¶
As mentioned above, because the keyword keep was used. The results are also stored in a dictionary. These data can
also be used for charting.
The reason for placing the results in a dictionary is to enable comparison of many scenarios, not just the first and the last.
with msolow.set_smpl(1,30):
msolow.keep_plot('income con*' );
msolow.modeldash('INCOME',jupyter=1)
apprun
Dash app running on http://127.0.0.1:5001/
Dash_graph(mmodel=<
Model name : Solow model
Model structure : Simultaneous
Number of variables : 11
Number of exogeneous variables : 5
Number of endogeneous variables : 6
>, pre_var='INCOME', filter=0, up=1, down=0, time_att=False, attshow=False, all=False, dashport=5001, debug=False, jupyter=1, show_trigger=False, inline=False, lag=False, threshold=0.5, growthshow=False)
3.8. More advanced topics¶
3.8.1. The logical stucture¶
Now the model has been analyzed, and the structure can be displayed.
You will find more on the logical structure here
3.8.1.1. Model structure¶
msolow.drawmodel( title="Model Structure", png=latex,size=(15,15))
3.8.1.2. Adjacency matrix¶
Another way to illustrate the dependency graph is an adjacency matrix.
msolow.plotadjacency();
The variables [‘INVESTMENT’, ‘CONSUMPTION’, ‘CAPITAL’, ‘INCOME’] in the red area are the core of the model and has to be solved as a system.
LABOR is the prolog and can be calculated before the core is solved. While CAPITAL_INTENSE is the epilog which can be calculated after the core is solved.
Many models comform to this pattern. And for solving purpose a model is divided into a prolog, core and an epilog. Even if the core is actually consistent of several strong components.
3.8.2. The python function used to solve the model¶
In order to solve the model Modelflow will generate a python function which implements the model. The user will hopeful newer have to relate to the generated python code. The point of modelflow is, that the user has to relate to the specification of the business logic, not the implementation in code
print(msolow.make_los_text)
def make_los(funks=[],errorfunk=None):
import time
import tqdm
from numba import jit
from modeluserfunk import jit, recode
from modelBLfunk import array, classfunk, clognorm, exp, gamma, inspect, jit, lifetime_credit_loss, log, logit, logit_inverse, lognorm, matrix, mv_opt, mv_opt_prop, norm, normcdf, qgamma, sqrt, sum_excel, transpose
def prolog0(values,outvalues,row,alfa=1.0):
try :
pass
values[row,0]=values[row-1,0]+(values[row,5]*values[row-1,0])
except :
errorfunk(values,sys.exc_info()[2].tb_lineno,overhead=9,overeq=0)
raise
return
def prolog(values,outvalues,row,alfa=1.0):
prolog0(values,outvalues,row,alfa=alfa)
return
def core0(values,outvalues,row,alfa=1.0):
try :
pass
values[row,9]=values[row,3]*values[row,1]**values[row,2]*values[row,0]**(1-values[row,2])
values[row,10]=(1-values[row,6])*values[row,9]
values[row,7]=values[row,9]-values[row,10]
values[row,1]=values[row-1,1]+(values[row,7]-values[row,4]*values[row-1,1])
except :
errorfunk(values,sys.exc_info()[2].tb_lineno,overhead=20,overeq=1)
raise
return
def core(values,outvalues,row,alfa=1.0):
core0(values,outvalues,row,alfa=alfa)
return
def epilog0(values,outvalues,row,alfa=1.0):
try :
pass
values[row,8]=values[row,1]/values[row,0]
except :
errorfunk(values,sys.exc_info()[2].tb_lineno,overhead=34,overeq=5)
raise
return
def epilog(values,outvalues,row,alfa=1.0):
epilog0(values,outvalues,row,alfa=alfa)
return
return prolog,core,epilog